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_ ^rhe  diffusion  model  of  reaction 
rates),  originally  due  to  H.  Kramers, 
is  rederived  and  extended.  The 
derivation-) follows  the  work  of  II' in 
and  Khasminskii  and^s  based  on  a 
clear  physical  picture  of  the  mole- 
cular events.  The  origin  of  the 
stochastic  forces  is  also  clearly 
treated.  Classical  mechanics  is 
used  throughout.  this  paper, 

jwe-  uses^the  Einstein-Snioluchowski 
approximation  and,  thus,  considers 
a diffusion  model  in  position  space 
only.  -We  non-dimensional the 
diffusion  equations* and  _pbta4-ft  a 
number  of^  singular  perturbation 
problems^  11  By  using  the  diffusion 
model,  one  can  treat  a number  of 
problems  involving  reaction  rate 
theory,  flc  3 npw  £prm  of 

transition  state  theory*.’  cal- 

culate reaction  rate  constants, 
transmission  coefficients  and 


the  lifetime  of  the  activated 
complex.!  Kramers 4 Result  is  t 
leading^ term  in  the  asymptotic 


expansion  of  the  rate  constant 
that  we  calculate.  We  show  how 
absorption  spectra  dan  also  be 
derived  by  use  of  the  diffusion 
model . 
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INTRODUCTION 


Almost  forty  years  ago,  H.A.  Kramers  (1)  introduced  a 
Brownian  motion  model  for  the  calculation  of  reaction  rate  constants. 
The  picture  involved  in  this  formulation  is  one  in  which  the  mole- 
cules undergo  a diffusion  process  in  reaction  space  while  moving 
in  physical  space.  Kramers'  theory  is  particularly  applicable  to 
reactions  in  solution,  and  reactions  of  relatively  large  molecule. 
Kramers  compared  his  theory  with  the  then  new  transition  state 
theory  (TST)  or  activated  rate  theory  associated  with  the  names 
Wigner,  Eyring,  and  Polayni.  The  Kramers  theory  has  been  relatively 
unnoticed  by  the  chemical  and  physical  community,  while  TST  and 
other,  more  complex  theories  have  developed  (2). 

As  a consequence,  there  are  excellent  classical,  semi- 
classical,  and  quantum  mechanical  methods  available  to  calculate 
rate  constants  for  the  reactions  of  small  molecules  in  the 
gaseous  phase.  On  the  other  hand,  when  large  molecules  or  polymers 
react  in  solution  a complete  theory  is  lacking.  The  Kramers 
diffusion  model  can  fill  this  gap (an  example  is  in  (3)).  In  this 
paper  and  the  following  one  (4),  we  extend  Kramers  theory  and  show 
how  the  diffusion  model  can  be  used  to  calculate  many  properties 
connected  with  reaction  phenomena,  not  just  reaction  rates. 

The  diffusion  approach  is  based  on  a Brownian  motion  model 
in  the  reaction  phase  space  of  the  molecule.  The  use  of  such 
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phase  space  distributions  was  initiated  by  Wigner  in  1932  (5) . 

By  the  use  of  Wigner' s technique  coefficients  in  the  differential 
equations  describing  the  evolution  of  the  phase  space  distribution 
function  can  be  given  approximate  quantum  interpretations  (6) . Our 
interest  here,  however,  is  the  solution  and  application  of  the 
diffusion  equations,  rather  than  derivation  of  the  equations. 
Consequently,  we  use  classical  mechanics  throughout.  The  extension 
to  Wigner  equivalent  formalism  is  straightforward  and  will  be 
considered  in  a later  paper  (7) . By  using  classical  mechanics, 
we  also  use  an  implicit  Born-Oppenheimer  assumption  about  the 
potential  energy  surface. 

Kramers  was  motivated  by  a desire  to  develop  an  alternative 
to  transition  state  theory  (1) . In  TST,  one  obtains  a reaction 

rate  constant  of  the  form: 

-Q/kRT 

k = k (f • f • ) e B (1.1) 

In  (1.1),  Q/k  T is  the  "activation  energy"  divided  by  Boltzmann's 
constant  times  temperature,  ( f - f * ) is  a frequency  factor  and  k 
is  the  transmission  coefficient.  Usually,  the  frequency  factor 
is  calculated  by  a quantum  mechanical  argument.  The  diffusion 
model  provides  a purely  classical  method  of  treating  the  frequency 
factor.  Hence,  we  obtain  CTST,  classical  transition  state  theory. 

CTST  and  TST  agree,  when  the  partition  functions  in  TST  are 
calculated  explicitly  (8). 


The  transmission  coefficient  is  usually  treated  as  an 


empirical  parameter.  TST  provides  no  method  of  calculating  k. 

The  diffusion  model  provides  a direct  way  to  calculate  k.  We 
will  show  that  TST  and  CTST  arise  as  special  cases,  of  the 
diffusion  model  (in  the  equilibrium  limit  and  vanishing  viscosity) . 

Previous  work  on  the  calculation  of  < was  done  by 
Hirschfelder , Wigner  and  Hulburt  (24)  using  quantum  mechanical 
techniques.  Here  a complementary,  stochastic  approach  is  given. 

Kramers  constructed  a stochastic  model  and  used  the  forward 
or  Fokker-Planck  diffusion  equation.  He  managed  to  construct 
solutions  of  the  Fokker-Planck  equation  in  certain  special  cases. 
The  deviation  between  TST  and  the  diffusion  theory  was  about  10-15 
percent  for  reactions  with  high  energy  barriers.  The  difference 
was  considerable  for  lower  barriers.  We  will  calculate  reaction 
rate  constants,  transmission  coefficients  and  lifetimes  of 
activated  complexes,  for  all  sizes  of  energy  barrier. 

In  section  2,  we  introduce  the  diffusion  model.  Our  analysis 
closely  follows  that  in  (9) . We  are  led  to  a system  of  stochastic 
differential  equations  in  the  reaction  phase  space.  Under  a 
second  limit,  the  Eins tein-Smoluchowski  (ES)  approximation,  the 
• structure  of  the  stochastic  differential  equation  is  simplified 
considerably.  We  give  the  forward  and  backward  equations  corre- 
sponding to  both  systems  of  stochastic  differential  equations. 

In  section  3,  the  equations  derived  in  section  2 are  non-dimen- 
sionalized.  A singular  perturbation  problem  arises.  In  section 
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4,  classical  transition  state  theory  is  derived.  In  section  5, 
we  use  the  ES  equations  to  calculate  reaction  rates,  transmission 
coefficients,  and  lifetimes  of  the  activated  complex.  In  section 
6,  we  show  how  spectra  can  be  calculated  using  the  diffusion  model. 
Not  all  of  our  results  are  new;  however,  we  are  unifying  many 
old  results  and  questions  with  a simple  conceptual  framework. 
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SECTION  2 


DIFFUSION  MODEL  AND  EINSTEIN-SMOLUCHOWSKI  (ES)  APPROXIMATION 

In  this  section,  we  introduce  the  diffusion  model  and  derive 
the  diffusion  equations  of  interest.  Our  approach  follows  II' in 
and  Khasminski  (9) . The  physical  assumptions  are  clearly  delineated 
and  the  introduction  of  stochastic  effects  is  also  obvious.  Since 
the  " is  based  on  dynamics,  one  could  use  the  projector 

operator  approach  (10) . However,  when  the  projection  operator  is 
used,  the  stochastic  assumption  is  introduced  in  a hidden  fashion, 
usually,  by  assuming  that  (e.g.,(lU) 

<-A(t)elLot  B(t)>  = <A(t)><B(t)>  (2.1) 

if  t > xc,  for  some  cut  off  time  tc  . In (2.1)  Lq  is  the 
Liouville  operator,  and  brackets  indicate  ensemble  averages.  For 
the  very  complicated  systems  of  interest  here,  it  is  not  possible 
to  use  dynamics  completely.  Namely,  it  is  not  yet  possible  to 
derive  the  stochastic  properties  from  the  dynamics  alone.  Hence, 
it  is  reasonable  to  introduce  the  stochastic  assumptions  at  the 
beginning.  The  work  of  II 'in  and  Khasminskii  (9)  is  such  a model. 

Since  their  paper  has  gone  unnoticed  by  the  physical  and  chemical 
community,  we  repeat  part  of  their  analysis  here  (equations  2.4  to 
2.16). 
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THE  DIFFUSION  MODEL 


Let  x denote  the  generic  reaction  coordinate  of  interest 
in  a molecule  of  mass  m . If  no  collisions  with  bath  molecules, 
of  mass  y , occur  then  classical  mechanics  is  obeyed,  so  that* 


x = 


mv  = F (x)  = -VV(x) 


(2.2) 

(2.3) 


In  (2.3),  V(x)  is  the  potential  function.  In  figure  1,  we  sketch 
V(x)  for  dissociation  reactions  and  for  tautomerizations . 

Let  £(t)  be  a stochastic  process  that  counts  collisions 
of  the  large  molecule  with  the  bath  molecules  and  let  F^  be 
the  distribution  function  of  £ : 


Pr(C (t)  < N } = F^ (N)  (2.4) 

Let  {t^}  i = 1,  2,  ...  be  the  jump  points  of  the  process  £ (t) . 

We  represent  the  bath  molecules  by  a family  of  identically  distri- 
buted random  variables  Uk },  k = 1 , 2 , . . . with  distribution 
function 

R (y)  = Pr( £k  < y)  (2.5) 


For  the  purposes  of  conceptual  simplicity,  we  will  treat  x as 
a scalar.  The  results  of  this  paper  immediately  generalize  to 
the  vector  case;  by  replacing  integrals  by  line  integrals. 


Introduce  a stochastic  process  (x^  (t)  , v (t)  ) , which  except  at 

the  points  r.  , coincides  with  (2.2,  2.3).  We  assume  that  x (t) 
K y 

is  continuous  and  v (t)  is  continuous  at  the  right  at  t.  , 

y x 

with  the  jump  at  given  by 


v*?] 


(2,6) 


Hence,  we  are  assuming  elastic  collisions.  If  we  consider  more 
than  one  space  dimension,  then  (2.6)  must  be  modified  by  the 
introduction  of  a factor  taking  into  account  the  spatial  distri- 
bution of  the  collisions  (4) . Define  the  transition  probability 


P^(x,  vi,  t,  xx,  v1)dx1dv1  = P^(x,  v,  t,  dxx,  dv^ 


(2.7) 


Pr{x^  (t)  Kx^Xj.  + dx^  , + dv^lx^o)  = x,  v^(o)  = v} 


Now  let  x(x,  v,  t),  v(x,  v,  t)  be  solutions  of  (2.2,  2.3)  with 
initial  conditions  x(o)  = x,  v(o)  = v.  Consider  a integrable 
function  f (x,  y)  and  define 

u (x , y,  t)  = rv  ^ (Xy (t) , vy(t)j 

= PrtTj^  > t}Ex^v[f  (xy  (t)  ,v^  (t) ) | tx  > t]  (2.8) 

t 00 

+ f dPr{ < S}ydR(z)Ex^v  [f  (xy,v^)  |t1=S,E1=z] 
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For  the  special  case  that  £ (t)  is  a Poisson  process  (the  case 
treated  by  II 'in  and  Khasminskii)  we  obtair. 


u(x,y, t)  = e~atf (x (x,v,t) ,v(x,v,t) ) 


t 


+f  ae  aSdS  J dR ( z ) u{x (x,v,S) ,v(x,v,S)  + v(z-v) 
*/0  — 


In  (2.9),  a is  the  parameter  of  the  Poisson  process  and 
v - 2y/(m+y)  is  the  reduced  mass.  We  note  that  x(x,v,t)  and 
v(x,v,t)  satisfy 

F 

y.  - vy  - — y = 0 

1 1 x nr  v 


dy  dx  , dv 
dt  dt  yx  dt  yv 


where  y = x or  y = v.  Thus,  if  f(x,v),  F(x)  are  three  times 
differentiable  (2.9)  becomes: 


9u 

at 


viii  + F(x)  m 
v3x  f lx;  av 


+ a 


oo 

y*u(x,v+v (z-y) ,t)dR(z) ; with  u(x,y,o)  = f(x,y) 


(2.9) 
,t-S} 

(2.10) 

(2.11) 

(2.12) 


Equation  (2.12)  was  derived  by  Il'in  and  Khasminskii.  Now 
consider  the  limit  a » 1 (i.e.,  many  collisions  per  second). 

Then  we  set 


where  n , T have  the  interpretations  of  the  viscosity  of  the 
medium  and  the  absolute  temperature,  respectively  (where  Boltzmann's 


We  assume  that 


constant  is 


The  assumptions  in  (2.14)  are  satisfied,  for  example,  by  the 


Maxwell  Boltzmann  distribution,  for  which 


Il'in  and  Khasminskii  prove  that  as  y 0 (2.12)  converges  to 


Equation  (2.16)  is  exact  and  rigorous  (compare  with  (12)).  It 


v dt 


(2.17) 


dx  = 


mdv  = (F  (x)  - nv)dt  + >J 2kBTn  dW 


(2.18) 


where  W(t)  is  the  Wiener  process.  We  call  2.16  - 2.18  the 
Ornstein-Uhlenbeck  (OU)  equations.  we  have  obtained  them  by 
introducing  two  stochastic  assumptions:  1)  a Poisson  process 
for  the  collisions;  2)  the  Maxwell  Boltzmann  distribution  for 
the  bath  molecules.  No  other  assumptions  are  needed. 

We  now  consider  a second  limit  of  (2.16),  the  Einstein- 
Smoluchowski  or  ES  limit.  In  the  limit  that  m/n+O  (high  vis- 
cosity) with  F/n  non-zero,  (2.16)  becomes  an  equation  for 
u(x,t);  independent  of  v: 


3u (x , t)  = afu  F(x)  3u^ 

at  n 3x2  n ax 


(2.19) 


Equation  (2.19)  corresponds  to  the  stochastic  differential  equation 


dx  = 


% 

_ F(X) 


dt 


2kBT 


dW 


(2.20) 


-12  -1 

For  molecules  of  molecular  weight  50,  m/n  ~ 10  sec  (13)  and 

-13  -1 

for  polymers  m/n  ~ 10  sec  (10).  Hence,  if  the  reactions  of 

-12  -1 

interest  have  rates  that  are  much  greater  than  10  sec  , the 
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ES  limit  will  be  a good  approximation.  Physically,  this  will 
often  be  the  case. 


II' in  and  Khasminskii  construct  asymptotic  solutions  of 
(2.16)  and  (2.17)  for  e = m/n  small.  Here  we  shall  be  interested 
in  a different  type  of  scaling  and  will  construct  different  asymp- 
totic solutions. 

FOKKER-PLANCK  AND  EXPECTED  TIME  EQUATIONS 

Denote  the  right  hand  sides  of  (2.16,  19)  by  Lu  and 

LpU  respectively.  These  operators  have  formal  adjoints  L* 

ES 

and  L*  Let  p be  a function  such  that 

ES 


pLu  - uL*p  = divergence 


(2.21) 


or 


pLEgu  - uL*ESp  = divergence 


(2.22) 


Then  p will  satisfy,  at  least  weakly,  the  Fokker-Planck  equations: 


pfc(x,v,t)  = 


-vpx  + 


m 


vv 


-[MH 


= L*p 


in  the  OU  case,  and 


§ - k(^) 


= L*  P 
“ ESf 


In  (2.23)  a subscript  indicates  differentiation. 
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(2.23) 


(2.24) 


in  the  ES  approximation.  Equations  (2.23,  24)  are  the  Fokker- 

Planck  equations  and  hold  if  certain  boundary  conditions  are 

2 2 2 2 

met  (e.g.  ,p-*0  as  x + v u-»0  as  x + v -•<»).  No 

expansion  of  the  master  equation  is  needed  (12,13).  Since  (2.16,19) 
are  rigorous,  (2.23,24)  are  also  rigorous,  provided  that  the 
boundary  conditions  are  satisfied  (14,15,16). 

Now  consider  the  following  definitions 
t(x,v)  = / tu. (x,v,t)dt  (2.25) 

*f)  t 


t(x) 


(x,  t)dt 


(2.26) 


Then  it  is  easy  to  show  that 


-u  (x , V , 00 ) 

V” 

3 2t  , vat  . /f  (x)  - nv^ 

at 

(2.28) 

m2 

^ 9*  \ m / 

av 

-u  (x ,°°) 

V 

n 

a2t  f (x)  at 
n 9x 

(2.29) 

Equations  (2.28,29)  were  derived  by  Klein  (17)  and  Weiss  (18) 
by  different  means  and  solved  in  some  particular  cases.  None 
of  the  problems  solved  here  were  treated  by  Klein  or  Weiss  ( i 9 ) . 

INTERPRETATIONS  AND  REMARKS 

We  have  not  specified  f(x,v)  or  boundary  conditipns  for 
any  of  the  equations  derived  in  section  2.2.  By  a judicious 
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choice  of  f and  boundary  conditions,  we  can  solve  many  problems 
by  using  (2.16-2.25).  For  example,  with  reference  to  figure  1, 

A 

u(x,v,t)  could  be  the  probability  that  a particle  has  x(t)  it  X, 
conditioned  on  x(o)  = x,  v(o)  = v.  Then  p(x,v,t)  would  be  the 
density  for  the  particle: 

p(x,v,t)dxdv  = Pr{ x (t) £ (x, x+dx)  , v (t) € (v, v+dv) } 

^ /v 

Finally,  t(x,v)  would  be  the  expected  time  to  reach  x(t)  a x, 
conditioned  on  x(o)  = x,  v(o)  = v: 

/V  A 

t (x , v ) = E min{ t:x (t)  = x|x(o)  = x,  v(o)  = v> 

In  this  case,  appropriate  boundary  conditions  are  for  (2.16): 


u(x,v,t)  = 1,  lim  u(x,v,t) 

X->— 00 

A 

u(x,v,o)  = 0 unless  x = x 


for  (2.23): 


fj p(x,v,t)dxdv  = 1 lim  p(x,v,t)  = 


for  (2.28): 

A 

t (x,v) 


2 2 
x + v 


3t 


0 lim  f-(x,v) 


X"*’—  * 


Analogous  interpretations  hold  for  the  ES  limit  equations,  and 
will  be  discussed  in  later  sections. 
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(2.30) 


(2.31) 


(2.32) 


(2.33) 


(2.34) 


i 


SECTION  3 


NON-DIMENSIONALI7.ATION  - THE  SINGULAR  PERTURBATION  PROBLEM 


We  will  now  introduce  scaled  variables  and  derive  the  non- 
dimensional  versions  of  the  equations  in  the  previous  section. 
First  consider  the  OU  equations  (2.16,23,28).  Let  Q denote  an 

A 

energy  (for  example  the  height  of  the  barrier  from  xq  to  x in 
figure  1)  and  let 

v = v1 

m 

F = A20  F, 

m 


x = 


_ v/Qm 


~2 


t = -t* 

n 


w = J - w' 

Y 


(3.1) 


Then  the  OU  equations  become 


k~T  ^2. 


3"  = v.9i!  + f'  — + _J_  _ v.iH 

3 t ' V 5x'  + f 3v'  Q 3(v,j2  3v 


(3.2) 


3£ 

at' 


V i!r 


Q 3 (v ' ) 2 


- - 5V.t(F'  - 


(3.3) 


— u (x ' , V 1 , °°)  * 


V a2t’ 

Q 3 (v ' ) 


+ F'iil  + v'^1  - V^i 

r 9v'  av ' aw* 


3x ' 


3v 1 


(3.4) 
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non-dimensional  OU  equations  are 


Next 


Equations  (3.5-10)  are  singular  perturbation  problems.  W 


assume  that 


where  V(x)  is  a non-dimensional  potential  function.  The 


SECTION  4 


CLASSICAL  TRANSITION  STATE  THEORY  (CTST) 


In  this  section,  we  show  how  TST  fits  into  the  above  frame 


work.  Since  our  results  are  derived  without  recourse  to  partition 
functions,  Planck's  constant  never  appears.  Hence,  we  call  these 
results  purely  classical  transition  state  theory. 

The  major  assumption  is  that  equilibrium  prevails.  The 
equilibrium  solution  of  (3.6)  is 


where  c is  chosen  so  that 


p (x,v)dxdv 


Note  that  H 


+ V(x)  is  the  Hamiltonian  of  the  classical 


well  around  x 


going  towards 


constant  is  then 


where  K is  the  transmission  coefficient,  which  will  be  calculated 


The  flux  across  x is 


-x 


vp(x,v)dv 


A 

-V(x)/e 


(4.4) 


To  obtain  j , we  must  divide  this  flux  by  the  number  of  particles 
in  the  well  at  X0»N0  • Since  p(x,v)  also  gives  the  particle 
density,  the  number  Nq  is 


«o  ■// 


f f -^(v2/2  + V (x)  ) 
p(x,v)dv  = J J e e dxdv 

WELL 


(4.5) 


We  now  replace  V(x)  by  a Taylor  expansion 


V(x)  = V(xQ)  + h V" (xq) (x  - xQ) 


and  let  the  limits  in  (4.5)  tend  to  +<*>  . We  obtain 


(4.6) 


r r -7(v2/2  + V (x  ) + >5V"(x  )(x  - x )2) 

No  1 J i e dvdx 


C 7 k IT  e ~V(xo)/e  2*e 

f2iTe\/-.ri7 -r--r  e = — — ■ 


-V(xQ)/e 


N^"(x0) 


(4.7) 


(4.8) 


Hence,  we  obtain 


yv"(x  ) -V(x)/e  V (x  ) / e 

j a e e 


In  (4.9),  we  recognize  V(x)  - V(xQ)  as  "activation  energy" 


(4.9) 
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of  the  reaction.  If  we  set  V(x  ) = 0,  then  the  reaction  rate 

o 


is : 


(xo}  e-V(x)/e 


2tt 


It  can  be  verified  that  the  usual  TST  (8)  gives  this  result,  if 


the  partition  functions  are  evaluated  explicitly.  Often  '4/"  (xq)  /2rr 
is  identified  with  the  "frequency  of  vibration"  in  the  well  about 
x . 


(4.10) 


SECTION  5 


ES  EQUATIONS:  REACTION  RATES,  TRANSMISSION  COEFFICIENTS, 
LIFETIMES  OF  ACTIVATED  COMPLEXES* 


In  this  section,  we  show  how  the  ES  equation  (3.8—3.10) 
can  be  used  to  calculate  much  of  the  desired  information  about 
a chemical  reaction. 


RATE  CONSTANT  BY  MODIFIED  KRAMERS  METHOD 

In  this  section,  we  calculate  the  rate  at  which  particles 

A 

pass  from  the  well  at  xf  to  the  barrier  peak  x . Our  method 
follows  that  of  Miller  (21)  and  Ludwig (16);  it  is  more  accurate 
than  Kramers  approach.  In  addition,  Kramers  derivation  is  based 
on  difficult,  somewhat  obscure  physical  arguments. 

We  seek  a solution  of  (3.9)  of  the  form 


p(x,t) 


-A  t 


(5.1) 


Then  each  o (x)  satisfies 
n 


-A  o 

n n 


.2 

3 o „ 

; — ? atr(F(x)o  ) . 

9x2  9x  n 


(5.2) 


kThe  results  in  this  section  generalize  to  the  multi-dimensional  case 
immediately  if  the  integrals  are  replaced  by  multiple  integrals. 

This  generalization  is  possible  because  the  ES  approximation  yields 
a gradient  deterministic  system.  Also  see  (19). 
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We  shall  calculate  the  lowest  eigen  value,  Xq  , with  boundary 
conditions : 


aQ(x)  = 0 


SSV <V  = 0 


(5.3) 


The  first  boundary  condition  corresponds  to  absorption  of  particles 
A 

at  x . The  second  insures  a constant  number  of  particles  at  xq  . 

(Kramers  assumed  this  also,  in  a disguised  form.)  The  rate  at 

A 

which  particles  reach  x from  x is  then 


j(x)  = JJ t pfc  (x , t ) dxdt  = XQ, 


(5.4) 


if  °Q(x)  is  properly  normalized  (as  it  must  be). 

We  integrate  (5.2)  once  and  use  the  fact  that  o (x  ) = 0: 

x o 


"Xq/"  o 

x 

o 

V (x)/e 


0(y)dy  - G 55J  co(x)  + vxo0 


(5.5) 


Since  e is  an  integrating  factor  for  the  right  hand  side, 

(5.5)  can  be  rewritten  as: 


-w!"  •• 


(y)dydS  = 


V (x)/e 


aoe 


(5.6) 


Now  consider  the  integral  I(S)  , defined  by 


I(S)  = / a (y)dy 

xo 


(5.7) 
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In  the  vicinity  of  x 


we  expand  a(y)  as 


°o(y)  = °o(xo)  e V(Y)/G[1  + egx(y)  + e2g2(y)  + . ..] 


where  g,  (y)  is  the  k order  correction  to  the  equilibrium 

A 

distribution.  The  density  in  (5.8)  does  not  vanish  at  x , as 

A g 

it  must  to  satisfy  aQ(x)  = 0 (it  does  satisfy  o^(xq)  = * 

CO  A 

Let  0 (y)  be  a neutralizer:  0 is  a C function,  0(x)=O 

A 

and  0 (y)  =1  if  y is  far  from  x,  e.g.. 


0 (y)  = 


if  y = x 

. , A i 

if  y < x - e 


for  some  n , and  is  smooth  in  between  x - e and  x . Instead 
of  (5.8),  we  use 


oQ(y)  = 0o(xo)e  V^y^e  0(y)  + o(e) 


Then 


I(S)  = J oo(xo)0(y)e"V(y)/edy  + 0(e) 


If  the  integral  is  evaluated  by  Laplace's  method,  we  obtain 


I (S)  ~ o (x_) e 


-V (x  )/e 


Hi  L-n2.l\  + 0(/e)  . 


Using  the  result  (5.12)  in  the  integral  equation  (5.6),  we 

A 

obtain,  for  x = x 
A 

x - 


-bf 


„V(s)  /ej„  . „ w / ire  _ V(xo)/e 

ds  CTo(xo)Y  2V"  (xo)  6 


= ^o^o1' 


(5.13) 


since  aQ(x)  = 0 . Hence, 


A 

x 


X „ f eV(s)/Eds 

TT  £ X 


(5.14) 


Using  Laplace's  method  to  evaluate  the  integral  in  the  denominator 
yields 


/2V"  (x  ) V (x  )/e 


TTE 


O r-V(x)/E-\/2|v"  (x) 


TT  £ 


= |\/v"  (xo)  Iv"  (x)  I e-V(x)/ee  V<xo)/e 


(5.15) 


(5.16) 


when  comparing  (5.16)  and  (4.10),  we  see  that  the  second  "frequency 
A 

factor"  |V"(x)|  does  not  appear  in  (4.10).  Equation  (5.16) 
is  also  Kramers  result  (but  derived  in  a different  fashion) . 

The  reaction  rate  constant  is  then  given  by 

k = (5.17) 
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I 


nee** 


REACTION  RATE  CONSTANT  BY  EXPECTED  TIME  FORMULATION 

In  this  section,  we  present  an  alternative  formulation  for 
the  rate  constant.  It  has  the  advantage  that  we  avoid  having  to 
use  the  asymptotic  analysis  which  assumed  e = kgT/Q  is  small. 
Consequently,  the  technique  of  this  section  works  for  moderate 
or  large  e (low  barriers)  as  well  as  small  e (high  barriers) . 


~ a _ A _ 

t(x)  = E{t:x(t)=x,x(s)<xfs<t|x(o)=x 

A 

x(t)  eventually  crosses  x) 


(5.18) 


I 


A 

at  which  molecules  reach  x from  x as 

o 


A 


_1 

UxQ) 


(5.22) 


where , 


t(xo) 


NV(s)/e 


/ 


-V(y)/e 


dydS 


(5.23) 


To  see  that  (5.22)  is  equivalent  to  (5.16),  we  assume  e « 1 and 
use  Laplace's  method  twice.  We  obtain 


t(xo) 


* 

/ ie 

/ Tre 

S 2 | V" (x) | 

V2lv"  (X oH 

eV (x)/ee-V (xq) /e 


(5.24) 


so  that 


X 


o 


e-v(x)/eev(xo)/e 


V" (x) |V” (xQ) 


(5.25) 


in  agreement  with  (5.16).  On  the  other  hand,  the  result  (5.21) 
is  much  more  versatile  than  the  eigenfunction  calculation. 

First,  £ need  not  be  small,  so  that  (5.21)  can  be  used  to  des- 
cribe catalysis  reactions  (an  application  would  be  to  the  system 
considered  in  (6)).  Second,  we  can  allow  for  the  experimentally 
true  fact  that  when  the  system  is  prepared,  not  all  the  molecules 
are  exactly  at  xq  . Instead,  there  is  a distribution  of  molecules. 
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G (x) dx 


Pr {x (o) 6 (x , x+dx) } 


(5.26) 


Then,  instead  of  l/t(xQ)  we  should  define  the  rate  constant  as 


Aq  = l/<t (x) > 


where 

<t(x)> 


(x) G (x) dx 


(5.27) 


is  the  ensemble  average  of  t(x). 

The  expected  time  formalism  is  much  more  versatile  than 
the  Kramers-Miller-Ludwig  approach.  In  figures  2 and  3,  we 

compare  the  CTST,  modified  Kramers  (eigenvalue)  and  expected 
time  formulations  for  the  rate  constant.  The  potential  used  for 
the  calculations  was 


V (x) 


ax 


(5.28) 


4 3/2 

for  which  Q = - x 

As  e decreases,  the  Kramers  and  expected  time  formulations 
converge,  as  is  expected  from  the  asymptotic  analysis.  The  CTST 

provides  a reasonable  estimate  of  the  rate  constant,  which  is 
remarkable  in  light  of  the  assumptions  used  to  derive  CTST. 

In  figure  3,  we  compare  the  three  theoretical  forms  of  the 
rate  constant  with  Monte  Carlo  experiments  (500  trials) . 
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THE  TRANSMISSION  COEFFICIENT 

Up  to  this  point,  the  transmission  coefficient  k is 
unspecified.  We  now  will  provide  an  exact  definition  of  the 
transmission  coefficient  and  will  show  how  to  calculate  it.  Let 
u(x)  be  defined  as 

/v  A# 

u(x)  = Pr{x(t)  crosses  x=x^  before  x=xfi|x(o)=x}  (5.29) 


Then  u(x)  satisfies  a stationary  version  of  (3.8) 


n 9 ^u  , _ - . 3u 

0 = e ~2  + F(x)  9^ 


(5.30) 


Ox 


with  boundary  conditions 


u (x  ) = 0 

o 


U(Xj)  = 1 


;5.31) 


The  transmission  coefficient  is  defined  as 
A 

< = u(x) 

Namely,  k is  the  probability  that  an  activated  complex  becomes 
a product  before  returning  to  the  reactant  state. 

The  solution  of  (5.30)  is 
V ( s'  /e. 


(5.32) 


l 


ds 


u ( X ) = 


/ 


(5.33) 


eV  (s)/e  ds 
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The  structure  of  the  potential  energy  surface  around  x will 

completely  determine  the  transmission  coefficient.  If  the  potential 

A , A , 

is  symmetric  about  x and  IV" (x) | is  bounded  away  from  zero,  then 
Laplace's  method  yields  (figure  4a),  as  expected 

k ~ *5  + o ( e ) (5.34) 

On  the  other  hand,  it  is  possible  that  k is  close  to  zero  (figure 

4b)  or  much  close  to  1 (figure  4c) , depending  upon  the  shape  of 

the  potential  surface  (figure  4).  However,  since  detailed  know- 

A 

ledge  of  the  potential  surface  around  x is  only  needed, it  should 
be  possible  to  calculate  k for  many  cases  of  interest. 

The  diffusion  model  thus  provides  a way  to  calculate  the 
transmission  coefficient,  which  was  previously  (in  TST)  treated 
as  an  empirical  parameter.  In  table  1,  we  compare  the  theoretical 
transmission  coefficient  with  Monte  Carlo  experiments,  for  the 
potential  (5.28) . 

LIFETIME  OF  THE  ACTIVATED  COMPLEX 

Many  reactions  proceed  according  to  a mechanism  in  which 
a true  reactant  complex  is  formed,  e.g., 

F + C~H  . - C H . F - C-H..F  + H (5.35) 

2 4 2 4 2 3 
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In  such  a case,  the  potential  energy  surface  will  have  a double 


minimum  structure  as  shown  in  figure  5.  The  complex  is  formed 
A 

in  the  well  at  x . The  lifetime  is  the  time  that  the  complex 

A 

remains  in  the  well.  We  pick  two  points  a,b  with  a < x , 
A 

b > x+  (figure  5) . The  mean  lifetime  of  the  complex  will 


satisfy 


i a2t  L . . at 

-1  = e — 7 + F(x)  r— 

3x^  d 


t (a)  = t (b)  = 0 


In  (5.36),  t(x)  is  the  expected  time  to  hit  a or  b,  given 
that  x(o)  = x. 

The  solution  of  (5.36,37)  is 


t ( 


*>  - u 


x S 

V(S)/e|*  Q-V(y)/e, 


e " " I e ' ' dydS 
a a 


/ evlsl/£dS 


IJ 

/ 


V(S)/e 


dS 


I 


v(s  )/e  / e-v<y)/cdyds 


A definition  of  the  lifetime  of  the  activated  complex  is  then 


Tac  * t(x) 
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(5.36) 

(5.37) 


(5.38) 


(5.39) 


l 


On  the  other  hand,  when  experiments  are  performed,  activated 
complexes  are  produced  according  to  some  distribution  H(x). 

Thus,  a better  definition  of  the  lifetime  of  the  activated 
complex  is 

x = E (t (x) ) (5.40) 

ac  H 


A 

+ 

t (x)dll  (x) 
x 


(5.41) 


AN  APPLICATION:  BIPROTONIC  PHOTOTAUTOMERS ION 

The  phenomenon  of  biprotonic  phototautomer ism  is  discussed 
in  (22).  A double  potential  well  is  present  (figure  6)*. 

In  this  case,  x represents  the  length  of  the  H - N hydrogen 

bond  (figure  6b)  in  7-azaindole  or  the  H - O bond  in  the  formic 
acid  dimer.  The  UV  and  green  flouresences  in  figure  6a 
correspond  to  7-azaindole. 

The  well  at  xq  corresponds  to  the  7-azaindole  dimer,  at 
x^  to  the  tautomer.  The  excitation  of  the  molecule  in  the  ground 
state  will  produce  a distribution  <}>(x)dx  of  molecules  in  the 
first  excited  state.  We  assume  that  when  the  excited  molecule 
reaches  xq(x^),  it  flouresces  with  probability  pQ(p^)  or 
decays  radiationlessly  with  probability  1 - pQ(l  ~ p^). 

‘Another  problem  involving  a double  potential  well,  which  can  be 
treated  by  these  methods  is  discussed  by  D.  Chandler  in  J.  Chem. 
Phys.  68:2959(1978). 
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u(x)  = Pr {molecule  reaches  xq  before  x-jJ  starts  at  x} 
then  u(x)  satisfies  the  backward  equation 


u = eu  + F(x)ux 
xx 


U (x  ) = 1 u (x.  ) = 0 

o 1 


(5.42) 


(5.43) 


(5.44) 


Assume  that  the  events  {reaching  xq(Xj)]  and  {flourescing 
from  xq(x^)}  are  independent.  Also,  we  assume  that  all  molecules 
with  x i xq  (s  x^)  reach  xq(x^)  and  flouresce  with  probability 


po (pl} * 

Let 

P u (x) 

I(X)  = Pl(l  - u(x)) 


(5.45) 


Then  I(x)  will  represent  a conditional  ratio  of 

flourescence  intensities.  The  total  flourescence 

x 


I 


I (x)  <)>  (x)  dx  + 


<t>  (x) dx 


X 


o 


<J>  (x ) dx 


UV/green 
intensity  is 


(5.46) 


The  first  term  in  (5.41)  is  the  contribution 
of  excited  molecules  initially  in  [xq,  x^.i 
the  contribution  from  molecules  with  x s x 


to  the  flourescence 
The  second  term  is 


1 


o * 


or  x £ x 


\%rfi  W 


SECTION  6 

CALCULATION  OF  SPECTRA 

Very  often,  we  are  interested  in  the  shape  of  the  "absorption" 
spectrum  of  a bond,  I (co)  . Let  <J>(t)  be  the  correlation  function 
of  x (t) : 

4>  ( T ) = E{x(t  + r)x(t)}.  (6.1) 


p (x,  o)  = 6 (x  - x) 


/' 


(6.5) 


p(x,  t)dx  = 1,  lim  plx’t)  > 0 

X -*a> 


Following  Ludwing  (16),  we  seek  a solution  of  (6.5)  in  the  form 


p(x,  t)  = e 


-iKx,t)/G  £ z 
k=0  k 


(6.6) 


with  4>(Xr  t)  and  Z^(x,  t)  k = 0 , 1 , 2 , . . . to  be  determined.  In 
many  cases,  it  is  sufficeint  to  use  only  the  first  term,  which 
will  be  accurate  to  order  e: 


p(x,t)  * Z e^(x't>/G  , 


(6.7) 


i.e.,  ZQ  is  a "normalization”  factor.  After  derivatives  are 
evaluated  and  substituted  into  (3.9),  terms  are  collected  according 
to  powers  of  e . The  leading  term  is  0 (e~ 'i'/' L / c)  and  vanishes  if 


ij)t  + F ( x ) - 0 (6.8) 

Ludwig  (16)  has  shown  how  this  equation  can  be  solved  by  the 
method  of  characteristics.  We  will  not  repeat  his  argument  here. 

We  note  that  (6.8)  corresponds  to  a "Hamiltonian" 

H = F (x) p + p2  (6.9) 
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— — 


and  to  "rays"  (where  p=^) 


at  ■ p<x’  + 2p-  af 


-Fx(xip  H - 1 


Along  these  rays 


(6.10) 


dip  2 
dS  = ? 


The  rays  cover  the  phase  plane.  Hence,  by  following  a sufficient 
number  of  rays  from  x for  a time  t , it  is  possible  to 
construct  (see  figure  7): 

p-(x,  x)dx  = Pr{x(x)6(x,  x + dx)jx(o)  = x} 

The  conditional  correlation  function  is  them 


(6.11) 


(6.12) 


mmm 
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APPENDIX  A 


THE  ENERGY  METHOD 


In  this  appendix,  we  briefly  describe  an  alternative 
formulation  of  the  problem  based  on  the  position-energy  phase 
space  rather  than  the  position-velocity  phase  space.  A special 
case  of  this  approach  was  used  by  Kramers  (1)  and  Visscher  (23). 
Our  approach  follows  Stratonovich  (24).  Consider  equations  (2.17, 
2.18).  The  total  energy  of  the  system  is 

2 

E = ^mv  + V ( x ) 


Hence 


dE  = mvdv  + V dx 
x 


Since  F(x)  = -v  , (A-2)  and  (2.17),  (2.18)  imply  that 

dx  = ^2/m(E  - v(x))  dt 

f 3 I 4k  Tn 

dE  = 2r)J-(E-V(x)  ) dt  + \ — 5 (E-V(x))  dW 

The  Fokker-Planck  equation  corresponding  to  (A-3,  A-4)  is 


pt (t,x,E) 


(A-l ) 

(A-2) 

(A-3) 

(A-4) 

(A-5) 


A-l 


In  ord^r  to  obtain  the  results  in  (1)  and  (23) , we  assume  that  the 


const  ^ E - v(x)  V (x)  < E 


otherwise 


Then  a simple  averaging  (as  in  (23) , page  117-119)  yields  an 


equation  for  the  density  p(t,  E)  in  E-space  only: 


(El  -]  , 2kBTn 
(E)  P m 


[W 

I Y * (E)  P 


- V(x)  dx 


(E  - v (x ) ) _Jjdx 


where  R(E)  = {x:  V(x)  < E}  . Equation  (A-6)  is  the  equation  of 
diffusion  in  energy  space  used  in  (1)  and  (23). 


(A-6) 


(A-7 ) 


APPENDIX  B 


MULTIPLE  BARRIERS 


In  the  body  of  the  paper,  we  did  not  consider  the  possibility 
that  more  than  one  energy  barrier  must  be  crossed.  The  results, 
however,  generalize  immediately.  For  example,  we  consider  the 

A 

expected  time  result  (5.23),  where  V(x)  now  has  k maxima  x 
j = 1,  ...,  k and  minima  xq^  , j = 1,  ...,  k . Then  the 
asymptotic  evaluation  of  (5.23)  is 


t(xo) 


1 V;  / 

£ Y'  / TTE 

e j 1 2 | V"  (x j ) 


_S£_V  eV<X1)/VV‘X«j,/e 

2V" ( x • ) / 

03  / 


(B-l ) 


where  xq  - x^^  . The  other  asymptotic  results  are  also  replaced 
by  expressions  involving  sums  over  various  contributions. 


CAPTIONS 


Figure  1:  Typical  potentials  of  interest  in  this  paper:  a)  a 
Lennard -Jones  like  potential;  b)  a double-minimum 
potential  (e.g. , hydrogen  bonded  tautomers) . 

Figure  2:  A comparison  of  numerical  results  using  CTST,  the  eigen- 
value formulation  and  the  expected  time  formulation  for 
the  reaction  rate  constant. 

Figure  3:  Comparison  of  CTST,  the  eigenvalue  formulation  and 
expected  time  formulations  for  the  reaction  rate 
constant  with  Monte  Carlo  experiments. 

Figure  4:  Transmission  coefficient  in  the  ES  formulation.  A 
steeped  barrier,  for  which  + 0(e)  > b)  a flat 

barrier,  for  which  k is  close  to  0;  c)  an  almost 
discontinuous  barrier,  for  which  k is  close  to  1. 

Figure  5:  Potential  energy  surface  for  a reaction  proceeding  by 
complex  formulation. 

Figure  6:  Potential  energy  surface  for  biprotonic  phototautomer ism. 

Figure  7:  A schematic  illustration  of  how  the  ray  method  can  be 
used  to  numerically  construct_correlation  functions. 

Rays  emanate  from  the  point  x . The  intersection  of 
the  line  t = t and  the  rays  determine  those  points 
that  can  be  reached  from  x in  time.  The  density  on 
these  points  is  p -(x,  r)  and  is  calculated  along 
the  rays. 

Table  1:  A comparison  of  the  theoretical  transmission  coefficient 

with  Monte  Carlo  experiments. 
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FIGURE  6 Continued 


TABLE  1 


'2500  Monte  Carlo  simulations  were  performed. 


COMPARISON 

OF  THE  THEORETICAL 
WITH  MONTE  CARLO 

TRANSMISSION 
EXPERIMENTS . 

COEFFICIENT 

Q(in  units 

of 

k0T) 

X 

< (Theory) 

K (MC)  * 

3 

1.31 

.58 

.58 

13 

2.14 

.53 

.53 

23 

2.58 

.52 

.52 
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